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Abstract 

This  paper  addresses  the  development  of  a  free  energy  model  for  magnetostrictive  transducers 
operating  in  hysteretic  and  nonlinear  regimes.  Such  models  are  required  both  for  material  and  system 
characterization  and  for  model-based  control  design.  The  model  is  constructed  in  two  steps.  In 
the  first,  Helmholtz  and  Gibbs  free  energy  relations  are  constructed  for  homogeneous  materials  with 
constant  internal  fields.  In  the  second  step,  the  effects  of  material  nonhomogeneities  and  nonconstant 
effective  fields  are  incorporated  through  the  construction  of  appropriate  stochastic  distributions. 
Properties  of  the  model  are  illustrated  through  comparison  and  prediction  of  data  collected  from  a 
typical  Terfenol-D  transducer. 
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1  Introduction 


Magnetostrictive  transducers  are  being  considered  for  a  number  of  high  performance  industrial,  au¬ 
tomotive,  biomedical  and  Department  of  Defense  applications  due  to  their  capability  for  generating 
large  force  inputs  over  a  broad  frequency  range.  As  detailed  in  [4,  5,  14],  applications  utilizing  mag¬ 
netostrictive  transducers  include,  among  others,  active  vibration  and  noise  control,  micropositioning 
in  high  force  regimes,  medical  and  industrial  ultrasonics,  non-contact  torque  sensors,  and  tuned  vi¬ 
bration  absorbers.  As  the  field  of  magnetostrictive  materials  becomes  more  mature,  it  is  anticipated 
that  a  growing  number  of  actuator  and  sensor  applications  will  benefit  from  the  multifunctional 
characteristics  and  robust  operation  that  these  materials  can  provide. 

An  inherent  property  of  all  existing  magnetostrictive  materials,  however,  is  the  presence  of  hys¬ 
teresis  and  constitutive  nonlinearities  in  the  relation  between  input  fields  and  the  generated  magne¬ 
tization  and  strains  as  illustrated  in  Figure  1.  While  the  effects  of  hysteresis  and  nonlinearities  can 
be  reduced  through  the  choice  of  stoichiometry  [21]  and  feedback  mechanisms,  the  modeling  of  these 
effects  in  a  manner  compatible  with  subsequent  control  design  can  yield  significant  improvements  in 
performance. 

To  illustrate  issues  pertinent  to  model  development  and  control  design  for  a  typical  magnetostric¬ 
tive  transducer,  we  consider  the  device  depicted  in  Figure  2.  Input  to  the  system  is  provided  by  a 
current  I(t )  applied  to  the  solenoid  which  generates  a  field  H(t).  This  causes  magnetic  moments 
in  the  Terfenol-D  rod  to  rotate  which  produces  the  strains  and  stresses  output  by  the  device.  The 
prestress  mechanisms  serve  two  purposes;  they  further  align  magnetic  moments  to  improve  perfor¬ 
mance  and  they  maintain  the  Terfenol  rod  in  a  state  of  compression.  The  surrounding  permanent 
magnet  provides  a  bias  field  Hq  necessary  to  achieve  bidirectional  strains  and  can  also  be  designed  to 
optimize  flux  paths  through  the  rod.  The  hysteretic  data  plotted  in  Figure  1  is  typical  for  this  type 
of  Terfenol-D  transducer  and  hence  it  represents  criteria  that  must  be  accommodated  in  models  to 
be  employed  for  design,  analysis  and  control  development. 

A  number  of  techniques  have  been  considered  for  quantifying  the  hysteresis  and  constitutive 
nonlinearities  inherent  to  the  relation  between  input  field  H  and  the  generated  magnetization  M 
and  strain  e.  For  Terfenol-D  transducers,  these  include  Preisach  models  [15,  17]  and  domain  wall 
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Figure  1.  Hysteretic  data  measured  in  a  Terfenol-D  transducer  as  reported  in  [6]:  (a)  field- 
magnetization  relation,  and  (b)  Field-strain  relation. 
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Figure  2.  Cross  section  of  a  typical  Terfenol-D  magnetostrictive  transducer. 


models  [1,  6,  7]  based  on  the  magnetic  theory  of  Jiles  and  Atherton  [9,  11].  For  transducer  char¬ 
acterization  and  control  design,  crucial  requirements  for  the  models  include  automatic  closure  of 
biased  minor  loops,  flexibility  with  regard  to  temperature  and  frequency,  and  sufficient  efficiency  to 
permit  real-time  implementation.  The  Preisach  techniques  guarantee  minor  loop  closure  but  require 
a  large  number  of  nonphysical  parameters  for  accurate  biased  minor  loop  characterization.  Fur¬ 
thermore,  modifications  to  the  classical  Preisach  theory  are  required  to  accommodate  temperature 
dependence,  reversible  changes  in  the  magnetization  observed  at  low  drive  levels,  or  to  relax  con¬ 
gruency  requirements  to  ensure  that  minor  loops  remain  inside  major  loops  [8,  12].  This  limits  the 
feasibility  of  control  techniques  based  on  Preisach  models  to  a  limited  number  of  applications.  Con¬ 
versely,  the  physical  basis  for  domain  wall  models  makes  them  applicable  to  a  broader  performance 
space.  However,  although  techniques  have  been  developed  to  guarantee  closure  of  biased  minor  loops 
when  turning  points  are  known  a  priori  [2,  10],  the  closure  of  minor  loops  when  these  points  are 
not  know  a  priori,  as  will  be  the  case  in  closed  loop  feedback  design,  is  not  guaranteed  by  present 
domain  wall  models.  This  proves  a  serious  detriment  in  control  design. 

In  this  paper,  we  quantify  hysteresis  and  constitutive  nonlinearities  inherent  to  Terfenol-D  mag¬ 
netostrictive  transducers  through  the  development  of  a  free  energy  model.  In  the  first  step  of  the 
development,  we  derive  appropriate  relations  for  the  Helmholtz  and  Gibbs  free  energies  at  the  lattice 
or  domain  level.  This  yields  a  local  magnetization  model  for  homogeneous,  single  crystal  materi¬ 
als.  The  effects  of  material  nonhomogeneities  and  nonuniform  effective  fields  are  then  incorporated 
through  consideration  of  stochastic  distributions.  This  yields  a  macroscopic  model  which  quantifies 
the  hysteretic  relation  between  H  and  M.  A  more  general  free  energy  relation  that  includes  magne¬ 
toelastic  interactions  is  then  employed  to  develop  nonlinear  constitutive  relations  which  predict  the 
nonlinear  relation  between  fields  H  and  strains  e. 

The  model  quantifies  hysteresis  in  the  relations  between  H  and  M  and  H  and  e  for  temperature 
invariant  and  quasistatic  operating  conditions.  Although  certain  relaxation  mechanisms  are  included, 
the  model  presently  neglects  eddy  current  losses  and  hence  it  should  be  employed  for  low  frequency 
regimes.  Furthermore,  the  model  is  developed  for  transducer  configurations  in  which  prestress  lev¬ 
els  are  sufficiently  high  to  dominate  crystalline  anisotropies  and  thus  in  this  initial  development, 
anisotropy  energy  is  neglected.  Finally,  it  is  illustrated  that  the  model  ensures  the  closure  of  biased, 
nested  minor  loops  in  both  the  H-M  and  H-s  relations. 

The  nonlinear  magnetization  model  is  developed  in  Sections  2.2  and  2.3  and  the  full  constitutive 
relations  are  developed  in  Section  2.4.  The  accuracy  of  the  model  is  illustrated  in  Section  3  through 
comparison  and  prediction  of  data  collected  from  a  typical  Terfenol-D  transducer.  The  capability  of 
the  model  to  ensure  closure  of  biased  minor  loops  is  illustrated  through  a  numerical  example. 
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2  Model  Development 


To  motivate  the  development  of  a  hysteresis  model  for  Terfenol-D,  we  summarize  first  pertinent 
properties  of  its  crystalline  structure.  This  is  facilitated  by  vector  conventions  in  which  planes  are 
denoted  by  round  brackets  and  directional  indices  are  represented  by  square  brackets.  Hence  the  faces 
of  a  cube  are  denoted  by  (100),  (010),  (001),  (100) ,  (010),  (001)  whereas  the  vertices  are  specified  by 
[100],  [010],  [001],  [100],  [010],  [001] .  In  both  cases,  I,  indicates  a  negative  direction.  Finally,  pointed 
brackets  are  used  to  summarize  an  entire  set  of  indices. 

In  present  manufacturing  processes,  Terfenol-D  crystals  are  grown  in  dendrite  sheets  oriented  in 
the  [112]  direction  as  depicted  in  Figure  3.  At  room  temperature,  the  easy  axes  lie  approximately 
in  the  (111)  set  of  directions  and  the  greatest  strains  occur  when  the  magnetization  M  rotates  from 
[111]  to  [111]. 

The  changes  in  magnetism  that  result  from  an  applied  field  H  are  primarily  due  to  two  mech¬ 
anisms:  the  rotation  of  moments  and  the  movement  of  domain  walls.  To  illustrate,  consider  a 
demagnetized  specimen  which  is  subjected  to  a  magnetic  field  oriented  in  the  [112]  direction  as 
shown  in  Figure  4.  At  low  field  levels,  the  change  in  magnetization  is  due  primarily  to  domain  wall 
motion  so  that  favorably  oriented  domains  are  enlarged.  As  the  field  is  increased,  moments  rotate 
to  orient  with  the  easy  [111]  axis.  This  produces  a  burst  region  in  the  H-M  or  H-e  curve  in  which 
small  changes  in  field  produce  large  changes  in  magnetization  or  strain.  In  the  final  stage  depicted  in 
Figure  4d,  the  material  acts  as  a  single  domain  as  moments  rotate  coherently  from  the  easy  axis  into 
the  direction  of  the  applied  field.  This  produces  the  saturation  behavior  exhibited  by  the  material. 

Strains  are  generated  by  the  material  when  moments  rotate  to  align  with  an  applied  field.  For 
general  configurations,  the  nragnetonrechanical  coupling  which  produces  these  strains  is  highly  com¬ 
plex  and  dependent  upon  factors  such  as  the  applied  stress  and  crystalline  anisotropies.  However,  in 
the  case  of  materials  in  which  a  prestress  perpendicular  to  the  moment  direction  is  sufficiently  large 
to  dominate  crystalline  anisotropy,  the  preferred  orientation  of  domains  is  shifted  from  the  original 
eight  (111)  magnetic  easy  axes  to  the  two  axes  [111]  and  [111]  perpendicular  to  the  [112]  direction. 
In  this  regime,  strains  are  due  primarily  to  moment  rotation  and  the  free  strain,  or  magnetostriction, 


Figure  3.  Orientation  of  Terfenol-D  crystals. 
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Figure  4.  Magnetization  process  in  the  (110)  plane  of  single  crystal  Terfenol-D  due  to  an  applied 
field  H  in  the  [112]  direction,  (a)  Demagnetized  state,  and  (b)  growth  of  domains  due  to  domain 
wall  motion,  (c)  Rotation  of  moments  to  the  easy  [111]  axis,  and  (d)  Rotation  of  moments  to  align 
with  the  applied  field. 


can  be  modeled  by  the  quadratic  relation 


A  = 


(1) 


where  Xs  and  Ms  respectively  denote  the  saturation  magnetostriction  and  magnetization.  As  dis¬ 
cussed  in  [1,  9],  the  prestress  levels  needed  to  optimize  transducer  performance  are  often  of  a  mag¬ 
nitude  such  that  stress  anisotropy  dominates  crystalline  anisotropy  and  the  relation  (1)  adequately 
models  the  strain  generated  by  the  material. 

To  model  the  hysteretic  behavior  of  the  material,  we  first  quantify  the  nonlinear  relation  between 
input  fields  H  and  magnetization  M.  We  then  develop  constitutive  relations  that  incorporate  the 
quadratic  behavior  (1)  to  provide  an  initial  model  for  the  strain  produced  by  the  transducer. 


2.1  Thermodynamic  Preliminaries 

To  provide  a  common  framework  for  constructing  magnetostatic  and  magnetoelastic  energy  relations, 
we  consider  general  properties  of  free  energy  relations.  The  free  energy  is  assumed  to  be  a  function 
of  temperature  T  and  an  order  parameter  e  which  is  respectively  taken  to  be  the  magnetization  M 
in  magnetostatic  relations  and  the  strain  e  in  magnetoelastic  relations.  We  let  (j)  denote  external 
fields  that  are  thermodynamically  conjugate  to  e.  For  the  order  parameters  M  and  £,  appropriate 
choices  for  the  external  field  4>  are  the  magnetic  field  H  and  stress  a.  Finally,  we  let  ip(e,  T)  denote 
a  general  Helmholtz  free  energy  relation. 

In  the  absence  of  applied  fields,  thermodynamic  equilibria  are  determined  by  minimizing  i/j  with 
respect  to  e  which,  under  the  assumption  of  differentiability,  yields  the  necessary  condition 

0(e’T)s^  =  o 

where  4>  denotes  the  energetic  response  of  the  system.  For  systems  subjected  to  an  external  field  (/>, 
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the  total  free  energy  is  taken  to  be 

^(e’T)  =  V>(e>T)  -  4>e  (2) 

which  yields  the  equilibrium  condition 

4>{e,T)  =  J  (3) 

The  relation  (3)  can  be  physically  interpreted  as  providing  conditions  under  which  the  order  pa¬ 
rameter  adjusts  to  balance  the  internal  energy  with  that  of  the  external  field.  As  will  be  illustrated 
in  subsequent  sections,  (3)  also  provides  the  constraints  used  to  specify  effective  fields  as  well  as 
magnetoelastic  constitutive  relations. 

2.2  Magnetic  Hysteresis  Model 

Motivated  by  the  prestressed  operating  regime  which  yields  the  [111]  and  [111]  easy  magnetic  axes, 
we  now  consider  the  energy  of  a  magnetic  moment  with  two  preferred  orientations.  This  includes  the 
characterization  of  the  internal  and  magnetostatic  energy.  Magnetoelastic  interactions  are  included 
in  Section  2.4  and  used  to  develop  elastic  constitutive  relations. 

We  first  formulate  a  Helmholtz  potential  J  by  considering  internal  energy  contributions  due  to 
moment  interactions.  These  interactions  are  assumed  to  be  isothermal  and  hence  temperatures  in  all 
subsequent  expressions  are  assumed  to  be  fixed  below  the  Curie  point,  which  implies  that  the  material 
under  study  is  in  its  ferromagnetic  phase.  This  condition,  when  combined  with  the  assumption  that 
stress  anisotropy  dominates  crystalline  anisotropy,  yields  material  behavior  consistent  with  a  double 
well  potential  of  the  form  illustrated  in  Figure  5.  As  summarized  in  [18],  statistical  mechanics 
analysis  determines  that  at  fixed  temperatures,  a  first-order  approximation  to  the  potential  exhibits 
quadratic  behavior  in  the  neighborhoods  of  all  three  equilibria.  Hence  we  employ  the  piecewise 
quadratic  definition 

(  ^{m  +  Mr)2  ,m<-Mj 

M )  =  |  lviM  -  Mr)2  ,  M  >  Mr  (4) 

[  -  Mr)  (fg  -  Mr)  ,  \M\  <  Mr 

for  the  Helmholtz  free  energy.  As  illustrated  in  Figure  5a,  Mj  and  Mr  respectively  denote  the 
inflection  point  and  magnetization  at  which  the  minimum  of  if}  occurs. 

From  (2),  the  Gibbs  energy  is  taken  to  be 

G  =  J-  HM  (5) 

which,  for  increasing  field  H ,  is  illustrated  in  Figure  5a.  Since  the  magnetostatic  energy  is  given 
by  £  =  no MH,  where  ho  is  the  magnetic  permeability,  a  second  choice  for  the  Gibbs  energy  is 
G  =  J  —  ho  HM.  However,  because  the  formulation  (5)  can  be  viewed  as  incorporating  ho  into 
the  two  energy  formulations  yield  equivalent  final  models,  and  we  employ  (5)  for  simplicity. 

For  a  homogeneous  material  with  effective  field  He  =  H,  where  H  denotes  the  applied  field,  the 
average  local  magnetization  is  given  by 

M  =  x+  {MJ)  +  x _  {MJ)  (6) 

where  and  x_  respectively  denote  the  fraction  of  moments  having  positive  and  negative  orienta¬ 
tions,  and  {MJ)  and  {MJ)  are  the  expected  values  of  the  resulting  magnetizations.  To  specify  {MJ), 
we  take 

p  OO 

{MJ)  =  /  Mn{G)dM  (7) 

J  Mo 
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V(M)=  G(0,M) 


Figure  5.  (a)  Helmholtz  energy  ip  and  Gibbs  energy  G  for  increasing  field  H  ( H-i  >  H\  >  0). 
(b)  Local  magnetization  M  as  a  function  of  H  for  a  homogeneous,  isotropic  material. 


where 

H(G)  =  Ce~GV/kT  (8) 

quantifies  the  probability  of  obtaining  the  energy  level  G.  k  is  Boltzmann’s  constant,  C  denotes  a 
constant  chosen  to  yield  a  probability  of  1  for  integration  over  all  possible  magnetization  values,  and 
Mo  is  the  critical  point  corresponding  to  the  unstable  equilibrium.  The  Boltzmann  energy  balance  is 
considered  over  a  lattice  volume  V  chosen  to  yield  relaxation  behavior  appropriate  for  the  material 
being  characterized.  Subsequent  evaluation  of  C  yields  the  average  magnetization  values 


We  note  that  when  implementing  the  model,  we  typically  replace  Mo  by  the  inflection  points  Mj 
and  —Mi,  respectively,  in  the  relations  for  (Af+)  and  {MJ).  This  simplifies  the  approximation  of  the 
integrals  and  can  be  motivated  by  observing  that  if  one  considers  the  forces  due  to  the  applied 
field,  maximum  restoring  forces  occur  at  Mj  and  —Mi  (e.g.,  see  pages  332-333  of  [3]).  Furthermore, 
for  materials  with  low  thermal  activation,  the  points  Mq  and  —Mi  coincide  in  the  limit  as  thermal 
activation  is  reduced  to  zero  for  positive  fields  while  M/  and  Mq  coincide  for  negative  fields  as 
illustrated  in  Figure  6. 

The  moment  fractions  x+  and  X-  are  quantified  by  the  evolution  equations 


x+  =  -p+-X. |_  +P-+X- 
X-  =  —p _ I- X -  +  p. I _ x+ 
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Figure  6.  (a)  Gibbs  energy  profile  with  a  high  level  of  thermal  activation  ( - )  in  the  Boltz¬ 

mann  probability  fx{G)  =  Ce~GV^kT .  (b)  Local  magnetization  M  given  by  equation  (6)  with  high 

thermal  activation  ( - )  and  limiting  magnetization  M  specified  by  (12)  in  the  absence  of  thermal 

activation  ( - ). 


where  m  is  the  mass  of  the  lattice  volume  V. 

Because  the  expression  (6)  for  the  local  magnetization  M  is  probabilistic  in  the  sense  that  the 
moment  rotations  are  determined  by  (8),  the  map  between  the  field  and  magnetization  exhibits 
both  hysteresis  and  nonlinear  transition  behavior  as  depicted  in  Figure  5b.  The  degree  to  which  the 
transition  is  mollified  is  dependent  on  the  ratio  between  GV  and  kT  in  the  Boltzmann  relation  (8); 
large  values  of  kT  model  regimes  in  which  thermal  activation  is  prominent.  This  in  turn  yields 
smoother  transitions  since,  for  a  fixed  field  level,  moments  have  a  higher  probability  of  achieving  the 
thermal  energy  required  to  overcome  energy  barriers. 

The  complex  behavior  of  M  can  be  simplified  in  two  aspects  to  facilitate  both  the  qualitative 
analysis  of  the  model  and  its  quantitative  implementation  for  regimes  in  which  thermal  activation  is 
negligible.  The  qualitative  behavior  can  be  ascertained  by  analyzing  the  equilibrium  behavior  after 
moments  have  switched.  In  this  regime,  the  equilibrium  condition  =  0  yields  H  =  J^.  Hence 
the  local  model  predicts  a  linear  relation  between  H  and  M  after  moment  switching.  Furthermore, 
from  =  Qpjji ,  it  follows  that  the  slope  of  the  hysteresis  kernel  in  this  linear  regime  is  K  The 
correspondence  between  the  critical  point  Mr  and  remanence  behavior  of  a  homogeneous  single 
crystal  follows  directly  from  the  zero  field  behavior  of  the  Gibbs  energy.  The  transition  point  occurs 
at  the  inflection  point  M j  since  this  is  the  point  of  maximum  restoring  force  [3].  In  concert,  these 
equilibrium  conditions  yield  criteria  which  can  be  employed  to  determine  initial  parameter  values  for 
optimization  routines  used  to  estimate  parameters  in  the  final  transducer  models. 
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For  regimes  in  which  operating  time  scales  make  thermal  activation  negligible,  asymptotic  re¬ 
lations  can  be  employed  to  simplify  the  magnetization  relations  since  jumps  in  this  case  can  be 
considered  to  occur  instantaneously  as  depicted  in  Figure  6.  In  this  regime,  M  =  Mmjn  where 
Mmin  is  obtained  through  the  solution  of  ^  =  0.  For  the  quadratic  Helmholtz  model  (4),  the  local 
magnetization  in  this  limiting  case  is  given  by 

(  [M(H;HC,£)}( 0)  ,  r(t)=0 

\M(H ;  Hc,£)](t)  =  |  f  -  Mr  ,  r(f)  /  0  and  H  (max  r(t))  =  -Hc  (12) 

{  f  +  Mr  ,  r(f)  /  0  and  if  (max  r(t))  =  Hc 

where  [M(H;  Hc,£)](0)  denotes  the  initial  moment  orientation  and  transition  points  are  specified  by 

r(t)  =  {t€  (0,  Tf]  I  H(t)  =  -Hc  or  H(t)  =  Hc}  .  (13) 

Employing  notation  used  for  Preisach  models  [20] ,  the  initial  moment  orientation  is  specified  by 

f  f  -Mr  ,  if (0)  <  -Hc 

[M(H-,Hc,0m  =  l  S  ,-Hc<H(0)<Hc  (14) 

{  f  +  Mr  ,  if  (0)  >  Hc  . 

The  dependence  of  M  on  the  local  coercive  field  Hc  =  t)(Mr  —  Mj)  is  explicitly  indicated  since 
parameter  i]  is  considered  to  be  distributed  for  the  bulk  model  developed  in  the  next  section. 

2.3  Nonhomogeneous  Materials  and  Effective  Fields 

The  local  magnetization  models  (6)  or  (12)  are  derived  under  the  assumption  that  the  lattice  and 
hence  domain  structure  in  the  materials  is  completely  homogeneous  so  that  free  energy  profiles 
for  different  regions  in  the  material  are  identical.  However,  this  typically  is  not  the  case  due  to 
material  defects,  nonuniformities  in  the  crystalline  structure,  and  polycrystallinity.  Furthermore, 
the  models  assume  that  the  effective  field  He  at  the  domain  level  is  exactly  the  applied  field  H 
and  hence  they  ignore  magnetic  interactions  or  Weiss  field  effects.  Through  the  incorporation  of 
appropriate  stochastic  distributions,  we  extend  the  local  model  to  obtain  a  bulk  magnetization  model 
for  nonhomogeneous  materials  with  nonconstant  effective  fields. 

To  accommodate  general  material  nonhomogeneities,  we  assume  that  materials  exhibit  a  distri¬ 
bution  of  free  energy  profiles  which  produce  variations  in  the  width  of  the  hysteresis  kernels  predicted 
by  (6)  or  (12)  as  depicted  in  Figure  7.  Nonhomogeneities  are  incorporated  by  assuming  that  the 
parameters  Mr,  Mj  or  Hc  =  tj(Mr  —  Mj)  are  manifestations  of  an  underlying  distribution  rather 
than  constants  as  assumed  in  Section  2.2.  For  this  initial  model,  we  consider  Hc  to  be  normally 
distributed  with  mean  Hc.  In  this  case,  the  total  magnetization  is  given  by 

M(H)  =  /  M(H;  Hc,  £)/ (Hc)dHc  (15) 

Jo 

with  the  density 

f(Hc )  =  C1e~(Hc-1Ic)2lb.  (16) 

Here,  Cj  and  b  are  positive  parameters  and  M  is  specified  by  (6)  or  (12).  We  note  that  in  (15),  the 
lower  limit  of  0  reflects  the  requirement  that  hysteresis  kernels  have  nonnegative  width.  Alternatively 
/  can  be  specified  as  a  log-normal  density  to  reflect  the  positivity  of  the  kernel  widths. 


Figure  7.  (a)  Free  energies  associated  with  nonuniform  moment  distributions,  (b)  Variations  in  the 
local  coercive  field  Hc  for  the  hysteresis  kernel  due  to  differing  free  energies. 


The  second  extension  addresses  the  incorporation  of  mean  field  effects  due  to  interdomain  cou¬ 
pling.  In  the  models  of  Jiles  and  Atherton  [9,  11]  and  subsequent  magnetostrictive  transducer  models 
[1,  6,  7],  this  coupling  was  modeled  by  an  effective  field  of  the  form 

He  =  H  +  aM  (17) 

where  a  is  a  constant  mean  field  parameter.  The  effective  field  (17)  is  analogous  to  the  Weiss 
mean  field  which  quantifies  interatomic  interactions.  Rather  than  assume  a  constant  interaction 
coefficient  a,  we  assume  that  the  effective  field  can  exhibit  variations  due  to  nonhomogeneities  in  the 
distribution  of  magnetic  moments.  To  incorporate  these  field  variations,  we  consider  the  effective 
field  to  be  normally  distributed  about  the  applied  field.  For  fixed  Hc ,  the  magnetization  in  this  case 
is  given  by 

/OO  _  0  - 

C2M(H;Hc,£)e-iH-H)  /bdU.  (18) 

-OO 

The  variations  in  the  effective  field  produce  domain  switching  in  advance  of  the  remanence  point  in 
accordance  with  observations  from  experimental  data. 

The  full  magnetization  model  for  nonhomogeneous  poly  crystalline  materials  with  variable  effec¬ 
tive  fields  then  becomes 

i*  OO  POO 

[M(H)\(t)  =  C  /  \M(H  +  H,Hc,£)}(t)e-n2/~be-{H'-Hc)2/bdHdHc  (19) 

with  M  again  given  by  (6)  or  (12).  We  note  that  for  implementation  purposes,  the  integrals  are 
truncated  using  high-order  Gaussian  quadrature  to  achieve  low-order  systems  which  facilitate  efficient 
implement  ation . 

Although  the  model  (19)  does  incorporate  certain  relaxation  mechanisms,  it  does  not  yet  incorpo¬ 
rate  eddy  current  losses  so  its  use  should  be  restricted  to  low  frequency  drive  regimes.  The  inclusion 
of  eddy  current  losses  and  variable  temperature  operating  regimes  to  accommodate  Ohmic  heating 
can  be  addressed  by  including  appropriate  terms  in  the  energy  relations. 
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2.4  Magnetoelastic  Constitutive  Relations 


The  Gibbs  relation  (5)  incorporates  the  internal  energy  and  magnetostatic  energy  exhibited  by 
isotropic  materials  at  the  domain  level.  However,  it  neglects  the  magnetoelastic  coupling  which 
provides  the  materials  with  magnetostrictive  capabilities.  Motivated  by  experimental  evidence  that 
suggests  a  quadratic  dependence  of  strain  on  magnetization,  we  consider  the  magnetostrictive  relation 
(1).  Although  this  relation  is  customarily  employed  in  magnetostrictive  transducer  design  due  to 
its  simplicity,  it  can  only  be  justified  theoretically  in  applications  where  material  stresses  are  such 
that  stress  anisotropy  overwhelms  crystalline  anisotropy.  In  this  case,  magnetoelastic  coupling  can 
be  incorporated  through  consideration  of  the  magnetoelastic  Helmholtz  free  energy  relation 

tpe(M,  e)  =  'ip(M)  +  ^Yms2  -  Ym1£M2  (20) 

and  corresponding  Gibbs  energy 

G(H,  M,  e)  =  V’(M)  +  -YMe2  -  Ym^eM2  -  HM  -  ae  (21) 

where  ij)  is  specified  by  (4).  Here  YM  denotes  the  Young’s  modulus  at  constant  magnetization  and 
7  is  a  magnetoelastic  coupling  coefficient. 

For  regimes  in  which  thermal  activation  is  significant,  the  local  magnetization  M  is  specified  by 
(6)  with  the  Gibbs  energy  relation  (21)  employed  in  the  integrals  (9)  and  (11).  For  the  limiting 
case  of  negligible  thermal  activation,  which  is  determined  through  solution  of  =  0,  the  local 
magnetization  is  given  by 


where  Hc  =  t/(Mr  —  Mj),  r  is  given  by  (13),  and 

—  j  ^ 


[M(H,e-,Hc,0m  ,  T(t)  =  0 

v~2Y-v1e  -  >  r^)  +  0  and  H (max. r(t))  =  -Hc 

i)-2r  -!7'  +  r,-ZYM1£  ’  r(^)  +  0  and  H (max r(t))  =  Hc 


h _ Mm 

■n-2YM~fs  r?-2y*7# 


,  H( 0)  <  -Hc 
,  -Hc  <  H( 0)  <  Hc 
_ 11 _ | _ mk’1 _  H({)\  >  H 

77-2 YMye  ^  ri-2YMje  ’  n  ' 

We  note  that  (22)  and  (23)  reduce  to  (12)  and  (14)  in  the  absence  of  strains. 

The  elastic  constitutive  relation  is  determined  from  the  equilibrium  condition 

dG 


de 


=  0 


which,  as  indicated  in  (3),  yields 


a  = 


d^e 


de 


M 


(22) 


(23) 


The  coupled  constitutive  relations  for  the  undamped  magnetostrictive  material  are  then  given  by 
<7  =  YMe  -  YMjM2 

roc  roc  _  (24) 

/  /  M(H  + H,£-,Hc,Oz  '  e~(  c~  c)  '  dHdHc 

1 0  J — oo 


M(H,  e)  =  C 


with  M  specified  by  (22)  in  the  absence  of  thermal  activation  or  (6)  with  G  given  by  (21)  when 
thermal  activation,  or  relaxation  mechanisms,  are  significant. 
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3  Model  Validation 


To  illustrate  properties  of  the  model,  we  present  two  examples.  In  the  first,  the  model  is  used  to 
characterize  and  predict  the  magnetization  and  strains  produced  by  a  Terfenol-D  transducer.  The 
second  illustrates  numerically  the  ability  of  the  model  to  close  biased  minor  loops. 

3.1  Experimental  Validation 

We  first  consider  the  characterization  of  magnetization  and  displacements  generated  by  a  typical 
Terfenol-D  transducer  having  the  configuration  depicted  in  Figure  2  in  response  to  a  quasistatic 
input  current  I(t)  to  the  solenoid.  Details  regarding  the  transducer  construction  and  manner  through 
which  the  experimental  data  was  collected  are  provided  in  [6].  We  consider  here  data  collected  at 
1  Hz  with  a  prestress  of  1  ksi  (6.9  MPa).  Two  input  levels  yield  the  moderate  and  high  drive  level 
magnetization  and  strain  data  plotted  in  Figures  1  and  8. 

The  parameters  Mr  =  3.7  x  104  A/m,  q  =  14,  Hc  =  300  A/m,  b  =  lx  108  A2/m2,  b  = 
8  x  108  A2/m2,  C  =  2.52  x  10~8  and  7  =  4.5  x  10~15  m2/A2  in  the  relations  (24)  were  estimated 
through  a  least  squares  fit  to  the  high  drive  level  data  producing  the  model  response  plotted  in 
Figure  8.  The  model  with  the  same  parameter  values  was  then  used  to  predict  the  moderate  drive 
level  behavior.  Because  we  considered  a  fixed  prestress  regime,  we  neglected  components  of  the 
magnetomechanical  coupling  in  the  constitutive  relations  (24).  Hence  the  bulk  magnetization  M 
was  computed  using  the  local  magnetization  relation  (22)  with  e  =  0  and  the  magnetostriction 
A  =  7 M2  was  used  to  model  the  measured  strains. 

Although  certain  discrepancies  exist  between  the  modeled  magnetization  and  strains  and  the 
experimental  measurements,  the  overall  behavior  predicted  by  the  model  is  sufficiently  accurate  for 
material  characterization  and  control  design.  The  contraction  in  the  low  field  magnetization  data  is 
typically  attributed  to  crystalline  anisotropies  inherent  to  Terfenol-D.  The  present  model  does  not 
incorporate  anisotropy  energy  and  hence  lacks  mechanisms  necessary  to  achieve  this  low  field  change 
in  concavity.  The  discrepancy  between  the  model  and  data  in  the  tips  of  the  strain  loop  indicates 
limitations  in  the  quadratic  magnetostriction  model  and  potential  mechanical  coupling  which  must 


Field  (A/m) 


(a) 


(b) 


Figure  8.  Experimental  data  ( - )  from  [6]  and  model  response  ( - ):  (a)  field-magnetization 

relation,  and  (b)  field-strain  relation. 
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Figure  9.  (a)  Magnetization  predicted  by  the  model,  (b)  strains  profile  predicted  by  the  model,  and 
(c)  input  field  H  to  the  model. 

be  incorporated  through  a  PDE  model  for  the  rod  based  on  the  constitutive  relations  (24).  The 
development  of  a  distributed  rod  model  and  the  computation  of  total  strains  and  displacements  can 
be  addressed  in  the  manner  detailed  in  [6,  7]. 

3.2  Biased  Minor  Loops 

To  illustrate  model  properties  under  asymmetric  minor  loop  operation,  the  field  plotted  in  Figure  9c 
was  provided  as  input  to  the  model  which  yielded  the  magnetization  and  strain  responses  plotted 
in  Figure  9a  and  9b.  The  parameters  in  the  model  were  taken  to  be  those  specified  in  the  example 
of  Section  3.1  which  were  obtained  through  a  least  squares  fit  to  high  drive  level  data.  Loops  1 
and  3  illustrate  biased  minor  loop  behavior  while  loop  2  illustrates  the  ability  of  the  model  to  enforce 
closure  of  multiply  nested  minor  loops.  Loops  4  and  5  illustrate  biased  behavior  leading  to  saturation. 
When  combined  with  the  experimental  results  presented  in  Section  3.1,  the  behavior  illustrated  here 
provides  the  model  with  substantial  flexibility  for  material  characterization  and  control  design. 
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4  Concluding  Remarks 


This  paper  presents  a  model  theory  to  quantify  hysteresis  and  constitutive  nonlinearities  in  mag- 
netostrictive  transducers  operating  in  moderate  to  high  drive  regimes.  By  combining  free  energy 
analysis  at  the  lattice  or  domain  level  with  stochastic  distributions  to  accommodate  material  and 
field  nonhomogeneities,  we  obtain  a  low-order  macroscopic  model  that  automatically  ensures  the 
closure  of  biased,  nested,  minor  loops.  In  its  present  form,  the  model  incorporates  relaxation  mecha¬ 
nisms  but  neglects  eddy  current  losses  so,  to  avoid  simulation  inaccuracies,  it  should  be  restricted  to 
low  frequencies  of  operation.  Errors  introduced  by  employing  the  model  at  high  frequencies  include 
phase  shifts  and  overestinration  of  magnetization  values.  Additionally,  although  stress  or  strain  ef¬ 
fects  have  been  included  in  the  free  energy  formulations,  the  accuracy  of  the  model  for  quantifying 
full  magnetoelastic  coupling  has  not  yet  been  established.  As  indicated  in  Section  3,  this  will  involve 
in  part  the  formulation  of  a  PDE  model  for  the  Terfenol-D  rod  and  subsequent  Galerkin  or  finite 
element  discretization  to  obtain  a  finite  dimensional  model.  The  model  also  neglects  self-heating 
in  the  transducer  as  well  as  crystalline  anisotropy.  A  number  of  present  transducer  designs  employ 
water  cooling  to  maintain  approximately  isothermal  conditions  and  require  prestress  levels  at  which 
stress  anisotropy  dominate  crystalline  anisotropy.  For  these  regimes,  the  model  accurately  quantifies 
the  low  frequency  dynamics  of  the  transducer. 

We  note  that  the  free  energy  framework  for  this  model  originated  in  the  context  of  SMA  [13,  16] 
and  was  recently  extended  to  piezoceramic  compounds  [19]  where  the  distributional  analysis  concern¬ 
ing  variable  coercive  and  effective  fields  was  added.  Hence  the  theory  provides  a  unified  framework 
for  modeling  hysteresis  in  a  number  of  representative  smart  material  systems.  Furthermore,  it  was 
illustrated  in  [20]  that  the  framework  can  be  employed  to  provide  an  energy  basis  for  Preisach  mod¬ 
els  with  the  following  important  difference;  temperature  and  relaxation  dependencies  are  manifested 
in  the  kernel  or  basis  for  this  model  whereas  they  enter  the  parameters,  or  measures,  in  Preisach 
expansions.  This  has  important  consequences  for  control  design  since  it  provides  the  potential  for 
eliminating  the  gain  scheduling  required  to  accommodate  changing  parameters  in  broadband  control 
regimes  or  systems  with  fluctuating  temperatures. 
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